%TESTPCHOL A test script of PCHOL.
%   A testing script that creates a small SPD matrix and verifies that the
%   results of PCHOL are correct.

% Test matrix
A = spconvert([...
    [ 8     1     6 ]; ...
    [ 3     5     7 ]; ...
    [ 4     9     2 ]; ...
    [ 9     9     1 ]; ...
    [ 7     1     6 ]; ...
    [ 3     8     7 ]; ...
    [ 5     9     2 ]; ...
    [ 9     3     3 ]; ...
    ]);
% Ensure that A is SPD
A = A*A' + speye(size(A,1));
n = size(A,1);

fprintf('Testing pchol. All errors below must be zero.\n');
for k = 1:n
    [Q,D,L] = pchol(A,k);
    Z       = sparse(k,n-k);
    E       = L*[[D Z];[Z' Q]]*L' - A; % Error
    fprintf('k = %2d  error = %.2e\n', k, full(max(abs(E(:)))));
end
